*************************************************************************************************
use  "${Gpath}data\\estimation_sample.dta", clear 
keep if car_int=="01" | car_int=="02"
gen category="A-PCI"

foreach i in $XVar0 $Xvar0_phy $Xvar $Xvar_phy  {
	cap gen `i'=0
}

keep share_dcbo_* cpfxcnesxyear sp_aa sp_mm sp_cnes overlap category ///
$XVar0 $Xvar0_phy $Xvar  m_30d team_size hospitalxdow hospitalxmonth car_int proc_rea sample_nonmissing morte category

append using  "${Gpath}data\\estimation_sample_surgery.dta", ///
keep(past_read share_dcbo_* cpfxcnesxyear sp_aa sp_mm sp_cnes overlap $Xvar  $XVar0 $Xvar0_phy $Xvar_phy  m_30d team_size hospitalxdow hospitalxmonth ///
sample_nonmissing proc_rea car_int morte dias_perm ///
share_age  ///
share_white share_black  share_other_race share_race_missing rais_merge category)

foreach i of varlist share_dcbo_* n_team* {
    replace `i'=0 if `i'==. 
}

* sample 
keep if sample_nonmissing==1 
keep if sp_aa>=2009 & sp_aa<=2020
drop if sp_aa==2020 & sp_mm==12 

gen emergency=real(car_int)>=2 

egen idd=group(category cpfxcnesxyear)
egen categoryxhospitalxdow=group(category hospitalxdow)
egen categoryxhospitalxmonth=group(category hospitalxmonth)

*****************************************************************************************************
* Elective vs. Emergency surgeries 
*****************************************************************************************************
estimates drop _all 
local h=0 

* emergency  
qui reghdfe  m_30d  share_dcbo_* $XVar0 $Xvar0_phy i.team_size overlap if emergency==1 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

* elective  
qui reghdfe  m_30d  share_dcbo_* $XVar0 $Xvar0_phy i.team_size overlap if emergency==0 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

esttab M* using "${Gpath}\\table\\elective_emergency.csv",   keep( *overlap ) replace star(* 0.1 ** 0.05 *** 0.01) se b(4) ///
staraux brackets   ///
compress   begin(";")  delimiter(";") end(";")  scalar(Mean N ) mtitles("Emergency" "Elective")


*** Robustness for elective admisions - with and without controls 
estimates drop _all 
local h=0 

qui reghdfe  m_30d  share_dcbo_* overlap if emergency==0 , ///
a(idd  ) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

qui reghdfe  m_30d  share_dcbo_*  overlap if emergency==0 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

qui reghdfe  m_30d  share_dcbo_* $XVar0   overlap if emergency==0 , ///
a(idd  ) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

qui reghdfe  m_30d  share_dcbo_*   $Xvar0_phy i.team_size  overlap if emergency==0 , ///
a(idd   ) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

qui reghdfe  m_30d  share_dcbo_* $XVar0 $Xvar0_phy i.team_size overlap if emergency==0 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth) cluster(sp_cnes) keepsing
qui sum  m_30d if e(sample)
scalar Mean=r(mean)
eststo M`h++', addscalars(Mean Mean )

esttab M* using "${Gpath}\\table\\elective_with_without_controls.csv",   keep( *overlap ) replace star(* 0.1 ** 0.05 *** 0.01) se b(4) ///
staraux brackets   ///
compress   begin(";")  delimiter(";") end(";")  scalar(Mean N ) 

*****************************************************************************************************
* Baseline: Other procedures and Multiple hypothesis testing 
*****************************************************************************************************
preserve 
keep if emergency==1 
levelsof category, local(G)
estimates drop _all 
local h=0 
foreach i in  `G'   {
local a=`h++'
qui reghdfe  m_30d  share_dcbo_* $XVar0 $Xvar0_phy i.team_size  if category=="`i'" , ///
a(cpfxcnesxyear hospitalxdow hospitalxmonth )  keepsing resid 
predict R`a'm_30d, resid 
label var R`a'm_30d "m_30d: `i'"
qui reghdfe  overlap  share_dcbo_* $XVar0 $Xvar0_phy i.team_size  if category=="`i'" & e(sample) , ///
a(cpfxcnesxyear hospitalxdow hospitalxmonth )  keepsing resid 
predict R`a'overlap, resid 
label var R`a'overlap "overlap: `i'"
}

rwolf2  (regress R0m_30d R0overlap ) (regress R1m_30d R1overlap ) ///
(reg R2m_30d R2overlap ) (reg R3m_30d R3overlap ) (reg R4m_30d R4overlap ) (reg R5m_30d R5overlap ) ///
(reg R6m_30d R6overlap ) (reg R7m_30d R7overlap ), ///
indepvars(R0overlap, R1overlap, R2overlap, R3overlap, R4overlap, R5overlap, R6overlap, R7overlap) reps(100) cluster(sp_cnes) 

forvalues i=0/7 {
global Pvalue`i'=e(rw_R`i'm_30d_R`i'overlap)	
}

* Final results 
global mtitles ""

levelsof category, local(G)
estimates drop _all 
local h=0 
foreach i in  `G' {
local a=`h++'
qui reghdfe  m_30d  share_dcbo_* $XVar0 $Xvar0_phy i.team_size overlap if category=="`i'" & emergency==1 , ///
a(cpfxcnesxyear hospitalxdow hospitalxmonth ) cluster(sp_cnes) keepsing

qui sum  m_30d if e(sample)
scalar Mean=r(mean)

scalar Pvalue=${Pvalue`a'}

eststo M`a', addscalars(Mean Mean Pvalue Pvalue)
global mtitles " ${mtitles}   `""`i'""'  "
}

esttab M* using "${Gpath}\\table\\other_procedures.csv",   keep(overlap) replace star(* 0.1 ** 0.05 *** 0.01) se b(4) ///
staraux brackets   ///
compress   begin(";")  delimiter(";") end(";")  scalar(Mean Pvalue N )  mtitles( ${mtitles}) 
restore 


*****************************************************************************************************
* Balance of other procedures 
*****************************************************************************************************
zscore $Xvar overlap 
zscore $Xvar_phy  share_age  ///
share_white share_black  share_other_race share_race_missing rais_merge 

preserve 
keep if emergency==1
*local cc=`f++'
* labels 
qui do "${Gpath}\setup\\labels.do"
* Patient characteristics 
global NAMES ""
global NAMES_Z ""
local h=0 
estimates drop _all 
foreach i in $Xvar  {
    dis as red "`i'"
local a=`h++'
qui reghdfe z_`i' share_dcbo_* z_overlap  , a(idd) cluster(sp_cnes )
estimates store M`a'
global NAMES "${NAMES}  M`a'="${`i'}"  "
}

coefplot (M*), xline(0) keep(z_overlap) aseq swapnames ///
coeflabels( ${NAMES} , notick  labgap(2)) eqlabels(, gap(5)) ///
graphregion(color(white) margin(zero )) ciopts(recast(rcap)) ///
ylabel(,labsize(small))  title("Panel A. Patient characteristics", pos(11) size(small)) msymbol(O) xlabel(-0.1(0.1)0.1) xlab(, nogrid) ///
xtitle("Standarized coefficients") xsize(3.2) ytitle("Dependent variable" " ")

graph export "${Gpath}\figures\\balance_patient_categoery`name'.pdf", replace 


* Physician characteristics 
global NAMES ""
global NAMES_Z ""
local h=0 
estimates drop _all 
foreach i in share_age ///
share_white share_black  share_other_race share_race_missing rais_merge $Xvar_phy    {
    dis as red "`i'"
	local a=`h++'
	qui reghdfe z_`i' share_dcbo_* z_overlap  , a(idd) cluster(sp_cnes )
	estimates store M`a'
	global NAMES "${NAMES}  M`a'="${`i'}"  "
}

coefplot (M*), xline(0) keep(z_overlap) aseq swapnames ///
coeflabels( ${NAMES} , notick  labgap(2)) ciopts(recast(rcap))  eqlabels(, gap(5)) ///
graphregion(color(white) margin(zero )) ylabel(,labsize(small))  title("Panel B. Physician characteristics", pos(11) size(small)) msymbol(O) xlabel(-0.5(0.25)0.5) xlab(, nogrid) ///
xtitle("Standarized coefficients") xsize(3.2) 

graph export "${Gpath}\figures\\balance_physician_category`name'.pdf", replace 
restore 


*****************************************************************************************************
* Variation in expertise overlap 
*****************************************************************************************************
qui reghdfe   overlap share_dcbo_*    if emergency==0 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth)  keepsing res(res_elective)

hist res_elective , frac  xtitle("Expertise overlap" "(residualized)") ///
title("Panel B. Elective admissions", pos(11) size(medium))  graphregion(color(white) margin(zero ))
graph export "${Gpath}\figures\\histogram_residual_elective.pdf", replace 

qui reghdfe overlap  share_dcbo_*   if emergency==1 , ///
a(idd  categoryxhospitalxdow  categoryxhospitalxmonth)  keepsing res(res_emergency)

hist res_emergency , frac  xtitle("Expertise overlap" "(residualized)") ///
title("Panel A. Emergency admissions", pos(11) size(medium))  graphregion(color(white) margin(zero ))
graph export "${Gpath}\figures\\histogram_residual_emergency.pdf", replace 